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ABSTRACT 


This report presents the complete analysis and formulas necessary 
to quantify the interference experienced by a generic satellite 
communications receiving station due to an interfering satel- 
lite. Both satellites, the desired as well as tlie interfering 
satellite, are considered to be in elliptical orbits. Formulas 
are developed for the satellite look angles and the satellite 
transmit angles generally related to the land mask of the 
receiving station site for both satellites. Formulas for 
considering Doppler effect due to the satellite motion as well 
as the Earth's rotation are developed. The effect of the 
interfering-satellite signal modulation and the Doppler effect 
on the power received are considered. The statistical formula- 
tion of the interference effect is presented in the form of a 
histogram of the interference to the desired signal power ratio. 
Finally, a computer program suitable for microcomputers such as 
IBM AT is provided with the flowchart, a sample run, results of 
the run, and the program code. 
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SECTION 1 
INTRODUCTION 


In the late seventies, the number of communications satellites 
in service increased at a rapid pace. Most of these satellites 
were geosynchronous. While interference is always experienced by 
the communications systems where they share frequencies of 
transmission or use channels adjacent to those of powerful 
transmitters, interference from geosynchronous satellites is 
becoming an increasingly important factor in establishing and 
operating a ground station. In the eighties the interference 
problem at the ground station has become worse because low Earth 
orbiter satellites have been added. These satellites are even 
more disruptive than those in geosynchronous orbits because of 
their low altitude above the Earth. The interference produced by 
all these satellites is generally intolerable to such a 
sensitive communications link as the Deep Space Network Systent 
at Goldstone, California. 

The interference experienced by the victim receiving station may 
be quantified in two different ways: one characterizes the 
interference in terms of statistics of the problem; the other 
simulates the interference phenomena with a computer. The 
statistical characterization of the interference problem suffers 
from the fact that it is quite difficult to statistically model 
the interfering signals and situations. This attempt introduces 
approximations and hand-waving in the modeling and, conse- 
quently, results in inaccurate as well as limited quantification 
of the interference effects. Simulation of the interference 
effects on a computer also is not trouble-free, because the 
computer must be run for considerable amounts of CPU time to 
generate the desired histogram. 
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There have been many efforts to quantify interference phenomenon, 
some quite good, but in general they are either too restrictive 
in their scope of applications (statistical-result programs) or 
too cumbersome to use (main-frame computer simulations). 
Computer simulations become more acceptable if they are 
transferred to a microcomputer where the user generally has 
exclusive use of the computer, and the CPU time does not enter 
the decision to run the simulation. 

Recently, microcomputers have become reasonably fast and their 
memory capacity has increased considerably; hence it has become 
feasible to write a simulation on a microcomputer that would 
normally be written on a mini or a main-frame computer. The 
following is an attempt to provide formulas and simulation for 
the histogram of interference experienced at a ground station. 
The histogram generates the mean and variance of interference 
along with such parameters as the cumulative distribution 
function (CDF), the probability density function (PDF), and the 
probability of the interference exceeding a given level. 

This is a complete simulation containing computation of orbital 
positions of the desired and interfering satellites using the 
classical orbital elements, calculation of the satellite antenna 
look angles for both satellites and look angles (elevation 
angles) at the desired-satellite ground- station antenna, and 
computation of Doppler effect due to the motions of the 
satellites and the Earth’s rotation; it also computes the 
interference-to-signal-power ratio, taking into account losses 
suffered by the links. After computing the interference-to- 
signal-power ratios, the program then computes the statistical 
quantities. The program is expected to be of general use to 
system designers and frequency managers in selecting the proper 
frequency under an interference scenario. 
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SECTION 2 

INTERFERENCE PARAMETERS 


2.1 Satellite Interference Problem Definition 

A satellite communications system generally has a satellite in 
an orbit and a ground station in communication with it. The 
downlinking from the satellite to the ground station is done at 
an assigned frequency and at a particular time slot in the orbit 
of the satellite. The ground-station antenna tracks the satel- 
lite and receives the telemetry or ranging signal. Tracking is 
necessary to maximize the gain of the antenna and consequently 
minimize the power required from the satellite for transmission. 

During the time period when the ground-station antenna, 
presumably a good-sized antenna dish with reasonable gain, is 
tracking the desired satellite, interference may result if an 
interfering satellite enters the antenna beamwidth. The magni- 
tude of interfering power received depends upon the following 
factors: (1) The frequency separation between the desired- 

satellite frequency and the interfering - satellite frequency, 

(2) the RF power radiated by the interfering-satellite antenna, 

(3) gain of the interfering-satellite antenna in the direction 


of the desired-satellite 

ground 

station. 

(4) gain 

of 

the 

desired-satellite ground 

station 

in 

the 

direction 

of 

the 

interfering satellite, and finally. 

(5) 

the losses 

of 

the 

desired-satellite ground 

system. 

it 

is well known that 

the 


magnitude of the interfering power alone does not determine the 
interference to the desired ground station, but it is the ratio 
of the interfering power to the desired signal power that 
quantifies the interference. Accordingly, we define the 
interference to signal power ratio (ISP) as: 


interference power at the desired satellite ground 

station antenna terminal 

ISP = desired satellite signal power at the ground station 
antenna terminals 
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In the following computer simulation, this ratio is computed in 
decibels and used as the system degradation measurement parame- 
ter. It is true that there are many computer programs for large 
machines and minicomputers that do something similar to what is 
proposed here, but there are only a handful, if any, programs 
available for microcomputers. 


2.2 Orbital Determination of Satellite Position 

Elliptical orbits are assumed for both the the interfering 
satellite and desired satellite. Using the classical theory of 
satellite orbits, one can compute the position of the satellites 
in a straightforward manner. No attempts will be made here to 
explain the theory behind orbital mechanics of a satellite; there 
are many excellent books and publications in the literature 
toward that end. We will use the already developed results and 
formulas for our purpose. 

The first step toward orbital position determination is to 
compute the position of the satellite in the orbital plane of 
the satellite. The x axis passes through the orbital perigee and 
the z axis is perpendicular to the orbital plane. The coordinates 
of the satellite in the orbital plane (x Q , y Q , z Q ) are 
obtained by using the classical theory of orbital mechanics of 
the satellites. The following steps are needed to generate the 
position of the satellite. 

Solve for the eccentric anomaly E of the satellite. Eccentric 
anomaly is obtained as follows. Locate the point where a line 
perpendicular to the semimajor axis of the orbit and passing 
through the position of the satellite intersects the circle that 
circumscribes the satellite orbit. A line from the center of the 
ellipse to this point makes an angle E (the eccentric anomaly) 
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with the semimajor axis. The eccentric anomaly of the satellite 
is obtained from the following equation: 


E - e sin (E) 


1 

a 


^398613 ■ 52 j 


1/2 



( 1 ) 


where a is the semimajor axis of the orbit, e is the eccentricity 
of the elliptical orbit (e = 0 is a circle), and t - t is the 
time of perigee passage. The radial distance of the satellite 
from the origin (the Earth's center), r Q , can be calculated 
from 


r Q = a (1 - e cos (E)) (2) 

This distance will become useful in computing the power received 
by the receiving station from the radiating satellite. The 
orbital plane coordinates can be computed by using the following 
equations : 

a( 1 - e 2 ) - r Q 


y Q = [(1 - e 2 ) ((a - r Q ) 2 + a 2 e 2 ) / e 2 ] 1/2 


It should be noted that the z coordinate is 0 because the 
satellite is in the orbital plane. 

The next step involves conversion of the orbital coordinates 
into the geocentric coordinate system. The plane of the orbit 
will also be defined within this system. In this system, the x 
axis always points towards the first point of Aries, and the z 
axis passes through the Earth's geographic north pole. This 
coordinate system translates as the Earth revolves around the 
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Sun, but it does not rotate. Figure 1 shows the orbit of the 
satellite around the Earth. Angular distance measured eastward 
in the equatorial plane is called the right ascension. The two 
points at which the orbit penetrates the equatorial plane are 
called nodes; the satellite moves upward through the equatorial 
plane at the ascending node and downward through the equatorial 
plane at the descending node. The right ascension of the ascend- 
ing node and the angle that the orbital plane makes with the 
equatorial plane i together locate the orbital plane with 
respect to the equatorial plane. 

To convert the orbital plane coordinates (x Q , y Q , z Q ) of 
the satellite into the geocentric coordinates (x^, y^, z^), 
we use the following equations: 

x^ = [cos (Cl) cos (to) - sin (Cl) cos (i) sin (to)] x Q \ 

- [cos (Q) sin (to) + cos (to) cos (i) sin (Q)] / 

y^ = [sin (Cl) cos (to) + cos (fl) cos (i) sin (to)] x Q \ (4) 
+ [-sin (H) sin (to) + cos (Cl) cos (i) cos (to)] y Q l 
z i = [sin (i) sin (to)] x Q + [sin (i) cos (to)] y Q J 

The next step is to compute the rotational coordinate system 
coordinates of the satellite. In this system, the x axis 
passes through the point at which the prime geographic 
meridian passes through the equator. The z axis coincides with 
the rotational Earth's axis. The transformation from the geocen- 
tric coordinates (x i , y^ z ^) to the rotational coordi- 
nates (x f , y , z r ) can be made in the following way: 

x r = x i cos (Cl e T e ) + y i sin (Cl Q T 0 ) 
y r - -Xj sin (!i e T e > + y t sin <S) e T e ) 
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Figure 1. Satellite orbital geometry. 


where 


A is the angular velocity of the rotating system, and 
© 

T 0 is the time elapsed since the x axis of the rota- 
tional system coincided with the x axis of the geo- 
centric system. 

T = 99.6909833 + 36000 . 7689*T^ 
e c 

+ 0.00038708*T C + 0.25068447*t (degrees) (6) 

where 

T = (Julian Day - 24150207)/36525 
and t is measured in minutes after midnight universal time (UT) . 


2.3 The Subsatellite Point 


This is the point at which a line joining the Earth's center to 
the satellite intersects the Earth's surface. This subsatellite 
point will be very useful in computing various angles and 
distances in the sections to follow. The subsatellite point will 
be characterized in terms of the north latitude, L , and west 

S 

longitude, l g . Using the rotational coordinates (x r , y r , 
z r ), we may compute the latitude and longitude mentioned above 
by using the following formulas (Ref. 1): 


subsatellite point 
latitude (degrees) 


L s = 90 - cos 1 [z r / r Q ] 


(7) 


i 
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subsatellite point _ , 

longitude (degrees) s 

-tan _1 [y r /x r ] 

180 + tan _1 [y r /|x r | ] 

90 + tan _1 [ |x r /y r | ] 
tan _1 [ |y r |/x r ] 

It should be realized that there will be one subsatellite point 
for the desired satellite and another one for the interfering 
satellite . 


x r _> 0 and y r _> 0 

x_ < 0 and y > 0 
r — ■* r — 

x r < 0 and y r <. 0 


( 8 ) 


x„ > 0 and y < 0 
r — ■* r — 


2.4 Visibility of the Satellites 

Satellite visibility determines when communications to and from 
the satellite occur. For interference to occur, the interfering 
satellite and the desired satellite must transmit downlink data 
(e.g., telemetry and ranging) at the same time. This implies 
that for interference to occur at the desired-satellite ground 
station, the ground station should be able to "see" the desired 
satellite plus the interfering satellite and at the same time 
the interfering-satellite ground station should "see" the 
interfering satellite. This is a necessary condition for 
interference, but obviously not sufficient in itself. More often 
than not, satellite transmissions are not scheduled to begin as 
soon as the satellite is above the horizon, but the transmission 
begins when the satellite is above a predetermined elevation 
angle. Primarily, this is because of the land mask around the 
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receiving station or some particular multipath condition around 
the receiving station. This elevation angle sets the limits of 
visibility of the satellite. The visibility limits of a given 
satellite may be different for the two ground stations. 

Figure 2 shows the interference geometry for the desired- 

satellite ground station. Let the desired-satellite ground- 

station latitude be L d and the corresponding west longitude be 

1^. Let L^ and 1^ be the latitude and west longitude of 

the interfering-satellite ground station. In the previous 

section, we computed the subsatellite point; using those 

formulas we may compute latitude and longitude pairs for the 

subsatellite points of the desired and interfering satellites. 

Let L sd and ^sd k e the latitude and longitude of the 

subsatellite point of the desired satellite and L . and 1 . 

si si 

the latitude and longitude of the subsatellite point of the 
interfering satellite. Using the geometry and notations for 
distances and angles used in Figure 2 we get: 

= cos -1 [cos (L d ) cos (L gi ) cos (l d ~l gi ) 

+ sin (L d ) sin (L gi >] 

n 0 = cos -1 [cos (L. ) cos (L . ) cos (l.-l .) 

/ J. 5 X 1 S 1 

+ sin (L. ) sin (L . ) ] 

X SI 

rig = cos -1 [cos (L d ) cos (L gd > cos (l d -l gd ) 

+ sin (L d ) sin (L gd >] 


( 10 ) 


[x 2 + v 2 + z 2 ] 1/2 

lx rd + y rd + z rd J 
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DESIRED 

SATELLITE 



Figure 2. Interference geometry. 
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( 11 ) 


r, = 


r„ = 


r„ = 


[r 01 + r e 2 - 2 r e r 01 cos 

lr 02 + r e 2 ‘ 2 r e C 02 003 <V) 1/2 


[ r 02 + r e 2 r e r 02 cos ^ n 2 *' ^ 


It should be noted that the d subscripts in the above formulas 
denote the desired satellite quantity while the i subscripts 
denote the interfering satellite quantity. The pertinent 
elevation angles can now be computed with ease. 


Elj - 

cos -1 

“ C 01 

sin 

(n 3 ))/r 1 ] \ 


ei 2 =■ 

cos -1 

t(r 02 

sin 

(n 1 ))/r 2 ] ( 

(12) 

E1 3 " 

cos -1 

[(r 02 

sin 

(n 2 ))/ r 3 l | 


E1 l- 

e1 2 ' 

and 

E1 3 

are known. 

they may be used in 


deciding when the downlink transmission from the satellite 
starts, i .e. , if 

El^ > a 1# then desired-satellite transmission begins 

El 3 > a^, then interfering-satellite transmission begins 

where aj^ and a 3 are the predetermined angles at which the 
desired-satellite and interfering-satellite transmissions begin. 
These angles basically depend upon the land mask of the respec- 
tive ground station. Thus the necessary condition for interfer- 
ence to occur at the desired-satellite ground station is that 
the following set of theoretical conditions results in "true": 


(El x > a x ) n (El 3 > a 3 ) fi (EL 2 > 0) 


(13) 
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The last term is necessary because the desired-satellite ground 
station should be seen by the interfering satellite to produce 
interference. If the answer to the above condition is "false," 
there is no interference possibility at the desired-satellite 
ground station. 


2 . 5 Antenna Angles 

To determine the power received at the desired-satellite ground- 
station antenna terminals from the desired satellite as well as 
the interfering power from the interfering satellite, we need 
the following of f-the-boresite angles of various antennas: angle 
© 3 off the boresite of the desired-satellite ground-station 
antenna to the desired satellite, angle © 2 off the boresite 
of the desired-satellite ground station antenna to the inter- 
fering satellite, angle © 4 off the boresite of the desired- 

satellite antenna to the desired satellite ground station, and 
finally, angle ©^^ off the boresite of the interfering- 
satellite antenna to the desired satellite ground station. We 
will assume that © 3 = © 4 = 0. This amounts to assuming 

that the antenna pointing systems on the desired satellite and 
the desired-satellite ground station are working perfectly and 
hence both these antennas are aligned properly. It will be 

assumed that the interfering-satellite antenna pointing system 
is working perfectly and, as a consequence, the interfering 
satellite antenna will be pointing to the interfering satellite 
ground station perfectly. Figure 2 shows the pertinent angles. 
The computation of the angles is done as follows: 

The distance x in Figure 2 can be computed by 
2 

x = r 0 [2(1 - cos (L^) cos (L^) cos (1^ - 1^) 

- sin (L i ) sin (L d ))] (14) 
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Using the law of cosines we get, 


©j = cos 1 [(r 2 + r 2 - x 2 )/2 r 2 r 3 ] (15) 

To compute © 2 we will need the distance r 4 in Figure 2. 
Since we know the coordinates of both the satellites, we can 
compute this distance quite easily. 

r 4 “ [<x rd - x ri )2 + <y r d - *ri> 2 + <z rd ' z ri )2) 

Using the law of cosines again we get, 

© 2 » cos -1 [ (r 2 + r 2 - r 2 )/2 r x r 2 ] (15) 

It should be noted that these antenna angles can be computed in 
many different ways, and some of these ways lend themselves to 
computer applications with ease and accuracy. 


2.6 Computation of Doppler Effect 

It will be assumed that the desired-satellite station will be 
tracking the Doppler effect introduced by its own (desired) 
satellite motion and hence the Doppler effect will not be 
considered for the desired satellite while computing the desired 
signal power. Doppler effect introduced by the interfering 
satellite motion will cause some variation in the interfering 
power received at the desired-satellite ground station and will 
be computed. The velocity component of the interfering satellite 
along the line r 2 in Figure 2 causes the Doppler effect at the 
desired-satellite ground station. It should be noted that this 
is not the only velocity causing the Doppler effect. The Earth's 
rotation also adds to the Doppler effect. Both these components 
are computed below. 
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The magnitude of the velocity vector of the interfering satellite 

V can be computed using the following formula: 

V = 398613.52 [(2/r Q2 ) - (1/aj)] (17) 

where aj is the semimajor axis of the interfering-satellite 
orbit . We need the component of V along r^ to compute the 
Doppler effect. Toward this end, we will need to compute the 
angle $ and the flight path angle 6 shown in Figure 2. Since 
we have the subsatellite points for both satellites, the compu- 
tation of these angles is not difficult. The distance y shown in 
Figure 2 can be computed from: 

y 2 = r 0 (l - 2 cos(n 1 )) 

Using the law of cosines (see Figure 2) we obtain 

* = cos- 1 [<<r 02 - c e ) 2 + r 2 2 - y 2 )/2 r 2 (r Q2 - r e ) 2 ] (18) 

The flight path angle for the interfering satellite would also 
be needed for computing the Doppler-ef fect-producing velocity 
component. The flight path angle 6 is shown in Figure 2. 

6 = tan -1 [e sin (4> 0 )//(l + 2 e cos (<t> 0 ) + e 2 )] (19) 

where <j> 0 is the angle between the semimajor axis and the 
line joining the focus of the orbit to the satellite. 

Thus the velocity component of the interfering satellite toward 
the desired-satellite ground station is given by 

V = V cos ( (ir/2) + 6 - $) (20) 

V is one Doppler-ef fect-producing velocity component; the 
c 

other component is due to the Earth’s rotation. These two 
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components add vectorially to produce the velocity that is 
responsible for the total Doppler effect. It should be noted that 
in this report, the Doppler due to the Earth’s rotational com- 
ponent will be computed only approximately. Assuming a spherical 
Earth with an average radius of 6378 km, one may compute the 
Earth's rotational speed at the equator to be about 0.4638 km/s. 

Figure 3 shows the geometry of Doppler effect due to the two 

velocity components. The vector along r 9 is the component V 

and the component along the line TP is the Earth's rotational 

component. These two components are at an angle a. The line TP 

is tangential to point T (the desired-satellite ground station) 

and the point P is directly east of T. The latitude of point P 

is the same as that of point T and the longitude of point P is 

slightly different than the longitude of point T (A is assumed 

to be small). The distance between S and P, d , may be 

sp J 

computed from 

d sp 2 = tx ri - r e cos (L d ) cos ((1 - A)l d )] 2 

+ ty ri ~ r e cos (L d ) sin ((1 - A)l d )] 2 (21) 

+ [ z ri ” r e sin (L d )]2 


The distance between T and P 



may be computed from 


d tp = [r e cos (L d } sin ( V 

- r 0 cos (L d ) sin ((1 - A)l d )] 2 
+ [r e sin (L d ) - r 0 sin ((1 - A) l d ) ] 2 


( 22 ) 


Using the law of triangles we compute the angle a 


a = cos' 1 [<4 + d* p - djj p )/(2 r 02 d tp >] 


(23) 
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The velocity component, V , which is due to the Earth's rota- 
tion and contributes to the Doppler effect, can be computed as 

V r = 0.4638 cos (a) (24) 

It should be noted that the velocity components V c and V r 
are along the line joining the desired-satellite ground-station 
antenna center to the interfering satellite and hence can be 
added algebraically. Thus the Doppler causing velocity, V^, is 
found from 


v d 

" V c + V r 




If the 

interfering 

satellite 

transmits 

at 

quency. 

then the frequency received at 

the 

ground 

station, f 

(which is 

modified 

by 

satellite motion), is 

calculated 

from: 



(25) 

f carrier fre- 
c 

desired-satellite 
the interfering- 


(1 + (V d /c)) 
fr [1 - (Vj/c) 2 ) ] 1/2 


(26) 


where c is the velocity of light. 


2.7 The Effect of Interfering Modulation 

The interfering power received by the desired-satellite ground 
station from the interfering satellite depends upon the modula- 
tion of the interfering signal and the front-end bandwidth of 
the receiver among many other things. Let the desired-satellite 
ground-station receiver front end have a bandwidth of BW Hz and 
the desired-satellite frequency be f d Hz. It is assumed that 
the receiver bandwidth, BW, is large enough for successful 
carrier tracking of the received doppler-corrupted desired 
signal. The receiver does not employ adaptive bandwidth 
optimization, i.e., once the bandwidth is seleted, it remains 
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constant. Figure 4 shows the baseband equivalent spectrum, and 
the shaded area shows the interfering power collected by the 
desired-satellite ground- station receiver. 


Interfering power 
collected by the 
receiver 



(PSD) df = FOF 


(27) 


where PSD is the baseband power spectral density of the modula- 
tion present on the interfering signal, delf is the frequency 
offset given by |f d - f r |, and E g is the energy in the 

pulse and P is the power in the pulse. For the NRZ pulse 
s 

modulation of the signal, E = P /data rate. 

s s 

It should be noted that the received frequency, f , used in 
the above calculation includes the Doppler effect. This implies 
that the shaded area will move according to the Doppler shift 
experienced by the incoming interfering signal. With all these 
parameters computed, we are ready to compute the interfering 
power received by the desired-satellite ground station and 
compare it with the desired signal power received. 


2.8 Power Received Calculations 

Suppose the desired satellite radiates P^ watts of RF power 
(this includes all the desired-satellite antenna inefficiencies). 
Using the previously defined parameters, the power input to the 
receiver, P dr , can be computed as 

P dr = P d L d g l (0) g 2 (0) tX 1 /(4irr 1 )] 2 (28) 

where signifies all the losses of the receiving ground station 
system (this may be frequency dependent). g 1 (0) is the gain of 
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POWER SPECTRAL DENSITY 



Figure 4. Interference modulation power spectrum. 
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the ground-station antenna and g 2 (0) is the 9 a ^ n of the 
,desired-satellite antenna. These antenna patterns are assumed to 
be symmetric in 4> direction. A similar equation for the 
interfering power received by the desired- satellite ground 
station may be written: 

P di = P i F0F L i g 3 (0 l ) g l (0 2 ) t x 2 /<4irr 2 )]2 (29) 


where is the interfering radiated RF power, signifies 
all the losses caused by the receiving ground station system to 
the interfering signal, g^O^ is the gain of the 
interfering-satellite antenna toward the desired-satellite 
ground-station antenna, g^©^ the gain of the 
desired-satellite ground-station antenna toward the interfering 
satellite, and FOF is the factor computed in equation (27) . 
Dividing equation (29) by the equation (28) we get a ratio ISR: 


I SR 


P i L i p2 r lf 
p d L a K r 2 J 


W gi< e 2> 

g i(°) g 2 (0) 


FOF 


(30) 


Normalizing the ISR by the factor (P^L^/P^L^) we obtain 
the previously defined interference measurement factor, ISP: 


ISP 



'MV 

g i(°) g 2 (o) 


FOF 


(30) 
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SECTION 3 

STATISTICAL CHARACTERIZATION 


In the theory described above we have assumed that the antenna 
pointing systems on the satellites as well as those on the ground 
stations function perfectly. This implies that the main lobe of 
the desired- or interfering-satellite antenna is always properly 
focused on the pertinent ground station as soon as the satellite 
becomes "visible" to the ground-station antenna. The interfering 
power received by the desired-satellite ground-station antenna 
depends upon the angles and 0 2< These two angles are 
functions of time and are constantly changing because of the 
change in the satellite orbital position. The interfering power 
also depends upon the Doppler effect on the interfering carrier, 
and this effect is also time dependent. Thus the received inter- 
fering power will change from orbit to orbit and from point to 
point in any particular orbit of these two satellites. Hence, 
even though the parameters of the interference scenario are con- 
stants, a histogram of the ratio ISP may be created giving rise 
to many statistical parameters of interest. There are many 
excellent books on the subject of histograms (Ref. 2); we will 
follow the procedure described below for constructing a 
histogram for ISP. 

Suppose we are interested in the range of ISP with the lower 
limit of ISPLL (dB) and the upper limit of ISPUL (dB) . Let N be 
the total number of sample ISP values gathered out of which 

N q be the number of ISP values such that ISP < ISPLL 

N 1 be the number of ISP values such that d Q < ISP < dj 
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N 2 be the number of ISP values such that 
N 3 be the number of ISP values such that 


d 1 < ISP < d 2 
d 2 < ISP < d 3 


be the number of ISP values such that ^k-1 < ISP < d k 

N k+1 1:56 t ^ ie num b er of ISP values such that ISPUL < ISP 
where 

k+1 

N = ^2 N i and d j = ISPLL + j DELISP, j = l,..,k (31) 

j=0 

DELISP = ^ ISPLL , = ISPLL> and = jspuL (32) 

The following figure depicts the bins produced by this process: 



N, 


N 


k-l 



N 


k+1 



ISPLL ISPUL 

Let be the middle point of the interval (dj_ 1 , dj); hence. 
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d. • , = ISPLL + DELISP ( j - 0.5 ) j = l,...,k+l (33) 

jmia 


1 k 

sample mean = m - jj N j d jmid 


= ISPLL + DELISP [0.5 + (( £ j N . )/N) ] 

j-1 D 


(34) 


1 2 
sample variance = jpy £ ^ N j^ d jmid “ ^ 


(35) 


sample probability = Pj = Prob C d. j < ISP < dj ] 


N, 

” j — 1 / • • • f ^+1 


(36) 


cumulative 
distribution = F. 
function 3 


k+1 

2 N i 


i=l 


N 


j = 1, . . . , k+1 


(37) 


sample 
probability 
density = p- 
function J 


N j 

N DELISP 


1 k+1 


(38) 


Using similar techniques we can also obtain the mean time 
between interference and average time the interference lasts. 
All these values are invaluable to the system designers and 
system analysts . 


It should be noted that the accuracy of the histogram depends 
upon the selection of DELISP and ISPLL as well as ISPUL. The 
information obtained above is sufficient to obtain some statis- 
tical parameters of interest such as the probability of inter- 
ference ISP exceeding a given ISP level . 
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SECTION 4 
COMPUTER PROGRAM 


The computer program is written specifically for a microcomputer 
such as an IBM AT. To make use of the program easier, it is 
written in IBM BASICA basic language. The program generally uses 
the theory and equations presented above; however, the theory 
does not discuss the complexity for computer implementation. The 
computer program harbors several changes to the equations 
described above to make them suitable for computer implementa- 
tion, always producing the same result as predicted by the 
theory. 

The basic flow of the computer simulation is to pass through a 
question and answer session for the user and gather all the 
parameters needed for the equations; this gives the user a fair 
amount of control over the execution of the program. All the 
computations are made and stored in double-precision arithmetic, 
but the printing is done in single precision. The inputs 
necessary for this program are listed below: 

Necessary inputs for the program: 

(1) Satellite orbital elements for the desired satellite. 

(a) Orbital eccentricity. 

(b) Semimajor axis. 

(c) Elapsed time since perigee. 

(d) Right ascension of the ascending node. 

(e) Argument of perigee. 

(f) Inclination of orbital plane to the equatorial 
plane . 

(g) Julian day of the first observation of the 
satellite. 
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(2) Satellite orbital elements for the interfering 
satellite . 

(a) Orbital eccentricity. 

(b) Semimajor axis. 

(c) Elapsed time since perigee. 

(d) Right ascension of the ascending node. 

(e) Argument of perigee. 

(f) Inclination of orbital plane to the equatorial 
plane. 

(g) Julian day of the first observation of the 
satellite. 

(3) Inputs for sampling the orbits. 

(a) Number of sample points per orbit. 

(b) Number of orbits to be sampled. 

(4) Desired-satellite ground-station parameters. 

(a) Gain of the ground-station antenna. 

(b) Ground-station receiver front-end bandwidth. 

(c) North latitude of the ground station. 

(d) West longitude of the ground station. 

(5) Interfering-satellite ground-station parameters. 

(a) North latitude of the ground station. 

(b) West longitude of the ground station. 

(6) Desired-satellite system parameters. 

(a) Frequency in use (MHz). 

(b) Data rate of the link (bits/s). 

(c) Gain of the antenna (dB). 

(d) Elevation angle at the ground station at which 
transmission begins (deg). 
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(7) Interfering-satellite system parameters. 

(a) Frequency in use (MHz). 

(b) Data rate of the link (bits/s). 

(c) Gain of the antenna (dB) . 

(d) Elevation angle at the ground station at which 
transmission begins (deg) . 

(8) Miscellaneous parameters. 

(a) Lower limit of range of interest of ISP (dB). 

(b) Upper limit of range of interest of ISP (dB). 

(c) Name of the file to store the histogram results. 

After these inputs are fed to the computer program, all the 
parameters necessary to evaluate IS/P at every sample point are 
computed and the following outputs are generated: 

Outputs of the program: 

(1) Average value of ISP (dB) . 

(2) Standard deviation of ISP (dB) . 

(3) Average duration of interference (minutes) . 

(4) Average duration free of interference (minutes). 

(5) Print the following outputs: 



(6) Store djmid' N j / P j » Pj / and Fj in the 
user-specified file, in the order specified above. 
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The program allows the user to enter any interfering-modulation 
spectrum. The user has to locate in the program the reminder 
statement, "Power spectral density of the interfering modula- 
tion" and starting from the next line the user may enter any 
interfering modulation. If nothing is entered then default is 
the spectrum of the NRZ pulse. Similarly, the user may enter any 
antenna pattern (symmetrical in 4> direction) for the desired- 
satellite ground-station antenna and the interfering-satellite 
antenna after the reminder statement, "Computing gain of the 
desired ground station" or "Computing gain of the interfering 
satellite" is found. 

It should be noted that there are 1000 bins made from the user- 
specified range of ISP, and the user does not have any control 
of this number. Hence the user should make sure that his storage 
device has enough space to store about 5000 double-precision 
numbers . 

The flow chart for the program is Figure 5, which shows the 
program’s steps to compute the ISP value at each sample point of 
the orbit. 
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Figure 5. Flowchart for the program. 
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PRINT, STORE DATA 


Figure 5 (contd) 
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SECTION 5 
AN EXAMPLE 


Let us suppose that we need to scope the interference environ- 
ment of a geosynchronous satellite communi cat ions system. The 
interference will be assumed to be generated by a low Earth 
orbiter satellite with a semimajor axis of orbit of 10000 km. We 
would like to compute the average value of the ISP (dB) , the 
standard deviation of ISP (dB), the maximum value of ISP (dB), 
the minimum value of ISP (dB), the average duration of interfer- 
ence experienced by the desired-satellite ground station in 
minutes, and the average duration in minutes that is free of 
interference . 

We would also like to compute the histogram of the interference 
(essentially the probability density of the interference) and 
the cumulative distribution function of the interference. The 
last curve would be helpful in finding such values as the 
probability of the interference exceeding a certain level. 

Figure 6 shows the inputs and Figures 7, 8, and 9 the outputs of 
the program for the example at hand. 
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SATELLITE ORBITAL ELEMENTS INPUTS FOR THE DESIRED SATELLITE mil 

ORBITAL ECCENTRICITY - O.OOOl 

SEMI -MAJOR AXIS (KM) - 42000.0 

ELAPSED TIME SICNCE PERIGEE (SECS) - 27916 

RIGHT ASCENT I ON OF ASCENDING NODE (DEG) - 70 

ARGUMENT OF PERIGEE (DEG) - 130 

ORBITAL PLANE INCLINATION (DEG) - 0 

THE JULIAN DAY OF FIRST OBSERVATION OF SATELLITE (DAYS) - 2443069.3 

SATELLITE ORBITAL ELEMENTS INPUTS FOR THE INTERFERING SATELLITE MM 

ORBITAL ECCENTRICITY • .0001 

SEMI -MAJOR AXIS (KM) » 10000.0 

ELAPSED TIME SICNCE PERIGEE (SECS) - 27916 

RIGHT ASCENTION OF ASCENDING NODE (DEG) - 50 

ARGUMENT OF PERIGEE (DEG) - 80 

ORBITAL PLANE INCLINATION (DEG) - 80 

THE JULIAN DAY OF FIRST OBSERVATION OF SATELLITE (DAYS) - 2443069.5 


ENTER THE NUMBER OF DESIRED SATELLITE ORBITS TO BE CONSIDERED AND THE NUMBER OF 
POINTS PER ORBIT - 5,300 


RECEIVING STATION PARAMETERS FOR THE DE6IRED SATELLITE inn 

ENTER GAIN OF THE ANTENNA (dB) * 40 
ENTER THE FRONT END BANDWIDTH - 2000.0 
ENTER THE NORTH LATITUDE (DEG) - 37 
ENTER THE WEST LONGITUDE (DEG) - 00 

RECEIVING STATION PARAMETERS FOR THE INTERFERING SATELLITE aaaaa 

ENTER THE NORTH LATITUDE (DEG) - 43 
ENTER THE WEST LONGITUDE (DEG) - 00 


DESIRED SATELLITE SYSTEM INPUTS aataiataii 

FREQUENCY IN USE (MHZ) - 1 
BIT RATE IN USE - 1000 


GAIN OF THE ANTENNA (DB) - 15 

ELEVATION ANGLE AT WHICH TRANSMISSION BEGINS (DEG) - S 


INTERFERING SATELLITE SYSTEM INPUTS aaaaaaaaaa 

FREQUENCY IN USE (MHZ) - 1 

BIT RATE IN USE - 1000 

GAIN OF THE ANTENNA (DB) - 15 

ELEVATION ANGLE AT WHICH TRANSMISSION BEGINS (DEG) - 5 


ENTER THE LOWER LIMIT OF RANGE OF INTEREST OF ISP (DB) - O 
ENTER THE UPPER LIMIT OF RANGE OF INTEREST OF ISP (DB) - 20 


ENTER THE FILE NAME TO STORE THE HISTPGRAM RESULTS • A BCD 


Figure 6. Input — output of the program. 
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ORBIT 

number 

1 

IS 

IN 

PROGRESS 

till! 

ORBIT 

NUMBER 

2 

IS 

IN 

PROGRESS 

1 t 1 1 t 

ORBIT 

NUMBER 

3 

IS 

IN 

PROGRESS 

1 I I 1 1 

ORBIT 

NUMBER 

4 

IS 

IN 

PROGRESS 

1 1 1 1 « 

ORBIT 

NUMBER 

5 

IS 

IN 

PROGRESS 

i t i i i 


% POINTS ABOVE THE SPECIFIED ISP UPPER LIMIT OK 20 <dB> - 
X POINTS BELOW THE SPECIFIED ISP LOWER LIMIT OF 0 (db) ■ 
MINIMUM VALUE OF ISP - -.630883 (db) 

MAXIMUM VALUE OF ISP - 18.90298 (dB) 

MEAN VALUE OF ISP - 11.15138 <dB) 

STANDARD DEVIATION OF ISP - 4.5681 <dB> 

AVERAGE DURATION OF INTERFERENCE - 16.57112 MINUTES 

AVERAGE DURATION FREE OF INTERFERENCE - 152.6125 MINUTES 


Figure 6 (contd) 


o 

1.075269 
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HISTOGRAM ENVELOPE 



Figure 7. Plot of the histogram envelope. 


I 
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All sn: 



0 5 10 15 20 

ISP(dB) 

Figure 8. Plot of sample probability density function. 
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ISP(dB) 

Figure 9. Plot of cumulative distribution function. 
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APPENDIX 

SYMBOLIC CODE OF THE PROGRAM 
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10 DEFDBL A-H,0-Z 
20 DEFINT I 

30 DIM IPDF(lOOl) ,POFY(1001) ,PDF(1001) ,DJMID(1001) ,CDF(1001) 

40 REM DATA 0.001181,42164.765,27916.3413,84.178,138.167,0.802, 
2443869.5 

50 REM DATA 0.1,10000.0,27916.3413,84.211,50,80.0,2443869.5 
60 REM DATA 5,10 

70 REM DATA 40,2000,37.229,80.438 

80 REM DATA 45,80 

90 REM DATA 1,1000,15,5 

100 REM DATA 1,1000,15,5 

110 REM DATA 5,20,ABCD 

120 FOR 1=1 TO 1001 : IPDF (I ) =0 : NEXT 

130 NDIV=100 : RE=6378 : ICONT=0 : ICOMT=0 : IORB=l : ICOUNT=0 :MSPECIAL=0 : 
NSPECIAL=0 : PI=3 . 1415926542# : CONVERT=PI/180 :TOTTIME56=0 : 
TOTTIME3 4=0 : PRINT 

14 0 TIME3=0 : TIME4=0 : TIME5=0 : TIME6=0 : IDB=0 : SUMDB=0 : PRINT "SATELLITE 
ORBITAL ELEMENTS INPUTS FOR THE DESIRED SATELLITE PRINT 

150 GOSUB 210 

160 EDSA=EEE : ADSA=AAA : TMTPDSA=TTT : CAPOMDSA=CCC : 0MDSA=000 : 

AINCLDS A=AI I I : AJDDSA=AJJJ 

170 PRINT "SATELLITE ORBITAL ELEMENTS INPUTS FOR THE INTERFERING 
SATELLITE PRINT 

180 GOSUB 210 

19 0 EISA=EEE : AISA=AAA : TMTPISA=TTT : CAPOMISA=CCC : 0MISA=000 : 
AINCLISA=AIII : AJDISA=AJJJ 
200 GOTO 300 
210 REM GOTO 280 

220 INPUT "ORBITAL ECCENTRICITY = ", EEE: INPUT"SEMI -MAJOR AXIS 
(KM) = " , AAA 

230 IF AAA<RE THEN INPUT "SEMI -MAJOR AXIS SHOULD BE GREATER THAN 
RADIUS OF THE EARTH, I.E., 6738 KM. PLEASE RE-ENTER : " , AAA 
240 INPUT "ELAPSED TIME SICNCE PERIGEE (SECS) = " , TTT 
250 INPUT"RIGHT ASCENTION OF ASCENDING NODE (DEG) = ",CCC: 
CCC=CCC* CONVERT : INPUT"ARGUMENT OF PERIGEE (DEG) = ",000: 
000=000 * CONVERT 

260 INPUT "ORBITAL PLANE INCLINATION (DEG) = ",AIII: 

AI I I=AI 1 1 * CONVERT 

270 INPUT"THE JULIAN DAY OF FIRST OBSERVATION OF SATELLITE (DAYS) 
= ",AJJJ: PRINT 

280 REM READ EEE , AAA, TTT , CCC, OOO, All I , AJJJ : REM PRINT EEE; AAA; TTT; 
CCC ; OOO ; All I ; AJJJ : CCC=CCC* CONVERT : 000=000* CONVERT : AIII=AIII* 
CONVERT 
290 RETURN 
300 REM GOTO 320 

310 PRINT :INPUT"ENTER THE NUMBER OF DESIRED SATELLITE ORBITS TO 
BE CONSIDERED AND THE NUMBER OF POINTS PER ORBIT = " , NUMORB , 
NPOINTS : PRINT 

320 REM READ NUMORB, NPOINTS : REM PRINT NUMORB ; NPOINTS 
330 TMTPO=TMTPISA: Cl=COS (CAPOMISA) : S1=SIN (CAPOMISA) :C11= 

COS (CAPOMDSA) : S11=SIN (CAPOMDSA) 

340 C2=COS (OMISA) : S2=SIN (OMISA) : C22=COS (OMDSA) : S22=SIN (OMDSA) 

350 C3=COS (AINCLISA) : S3=SIN (AINCLISA) : C33=COS (AINCLDSA) : 

S33=SIN (AINCLDSA) 

360 FAJDISA=AJDISA-INT (AJDISA) : FAJDDSA=AJDDSA-INT (AJDDSA) : 
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0RBPD=2*PI*ADSA*SQR (ADSA/398613 . 5) : TIME=0 : TIME6=0 : 
DELT=ORBPD/NPOINTS : SSTOP=NUMORB*ORBPD 
370 PRINT "RECEIVING STATION PARAMETERS FOR THE DESIRED 
SATELLITE PRINT 

380 REM GOTO 400 

390 INPUT "ENTER GAIN OF THE ANTENNA (dB) = " , GNDST : 

INPUT"ENTER THE FRONT END BANDWIDTH = " , BW 
400 REM READ GNDST, BW: REM PRINT GNDST ;BW 
410 GOSUB 470 
420 DSTLAT=DDD : DSTLON=EEE 

430 PRINT "RECEIVING STATION PARAMETERS FOR THE INTERFERING 
SATELLITE PRINT 

440 GOSUB 470 

450 AISTLAT=DDD:AISTLON=EEE 
460 GOTO 510 
470 REM GOTO 490 

480 INPUT" ENTER THE NORTH LATITUDE (DEG) = " , DDD: DDD=DDD*CONVERT: 
INPUT "ENTER THE WEST LONGITUDE (DEG) = " , EEE : EEE=EEE* CONVERT : 
PRINT 

490 REM READ DDD, EEE : PRINT DDD; EEE : DDD=DDD* CONVERT : EEE= 


EEE* CONVERT 
500 RETURN 

510 GNDST=10 A (GNDST/10) 

520 PRINT: PRINT"DESIRED SATELLITE SYSTEM INPUTS :::::::::: n : PRINT 
530 GOSUB 590 

540 FRDST=FFF*1000000 l : GNDSA=GGG : ANGDS A=AAA : DSTBR=BBB 

550 PRINT" INTERFERING SATELLITE SYSTEM INPUTS :::::::::: ": PRINT 

560 GOSUB 590 


570 FRIST=FFF*1000000 ! :GNISA=GGG:ANGISA=AAA: AISTBR=BBB 
580 GOTO 640 
590 REM GOTO 620 

600 INPUT" FREQUENCY IN USE (MHZ) = ", FFF: INPUT" BIT RATE IN 
USE = " , BBB 

610 INPUT "GAIN OF THE ANTENNA (dB) : " , GGG : IN PUT "ELEVATION ANGLE 
AT WHICH TRANSMION BEGINS (DEG) = ", AAA: AAA=AAA* CONVERT: PRINT 
620 REM READ FFF, BBB, GGG, AAA: PRINT FFF ; BBB ; GGG ; AAA : AAA= 

AAA* CONVERT 
630 RETURN 

640 GNDSA=10 A (GNDSA/10) : GNISA=10 A (GNISA/10) : GNIST=10 A (GGG/10) 

650 REM GOTO 690 

660 PRINT: INPUT "ENTER THE LOWER LIMIT OF RANGE OF INTEREST OF 
ISR (dB) = " , AISPLL 

670 INPUT "ENTER THE UPPER LIMIT OF RANGE OF INTEREST OF ISR 
(dB) = " , AISPUL: DELAISP= ( AISPUL-AISPLL) /1000 
680 PRINT: INPUT "ENTER THE FILE NAME TO STORE THE HISTOGRAM 
RESULTS = " , FILEA$ 

690 REM READ AISPLL, AISPUL, FILEA$ : PRINT AISPLL; AISPUL; FILEA$ : 
DELAISP=( AISPUL-AISPLL) /1000 

700 PRINT : XC=C1*C2-S1*S2 *C3 :XS=-C1*S2-S1*C2*C3:YC=S1*C2+C1*S2*C3 : 
YS=-S1*S2+C1*C2*C3 : ZC=S2*S3 : ZS=C2*S3 
710 XC11=C11*C22-S11*S22*C33:XS11=-C11*S22-S11*C22*C33: YC11=S11* 
C22+C11*S22*C33: YS11=-S11*S22+C11*C22*C33 :ZC11=S22*S33: 
ZS11=C22*S33 

720 PRINT: PRINT: PRINT: PRINT 
730 REM CONTINUE 
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740 ISEDSADST=0 : ISEISAIST=0 : ISEISADST=0 

750 DISDSTIST=RE*SQR ( 2 * ( 1-COS (AISTLAT) *COS ( DSTLAT) * 

COS ( AISTLON-DSTLON) -SIN (AISTLAT) *SIN( DSTLAT) ) ) 

760 IF (TMTPISA-TMTPO) > (SSTOP-DELT/600) THEN GOTO 1870 
770 IF (TIME> (IORB-1) *ORBPD/ 60) AND (TIME< ( IORB*ORBPD/60 ) ) 

THEN IORB=IORB+l: PRINT: PRINT"ORBIT NUMBER "IORB-1; 

••IS IN PROGRESS 

780 REM PREPARATION FOR SUBROUTINE COM1+COM2 FOR DESIRED 
SYSTEM :::::::::: 

790 AAS=ADSA : TMTPS=TMTPDSA : TIMES=TIME2 : EPSIS=EDSA : XCS=XC1 1 : 

XSS=XS11 : YCS=YC11 : YSS=YS11 : ZCS=ZC11 : ZSS=ZS11 : FAJDS=FAJDDSA 
800 AJDS=AJDDSA : DLATS=AISTLAT : DLONS=AISTLON : VLATS=DSTLAT : 
VLONS=DSTLON 
810 GOSUB 940 

820 DELTADSA= DELTAS : R0D=R0S : DISDSAIST=DISDS : DISDSADST=DISVS : 
ELISTDSA=ELEDS 

830 ELDSTDSA=ELEVS : GMISTDSA=GAMAD : GMDSTDSA=GAMAV : XR2=XR : YR2=YR : 
ZR2=ZR: DDDD=DDD2 

840 IF ELDSTDSA >= ANGDSA THEN ISEDSADST=1 
850 REM PREPARATION FOR SUBROUTINE COM1+COM2 FOR 
INTERFERING SYSTEM ::::::: 

860 AAS=AISA : TMTPS=TMTPISA: TIMES=TIME : EPSIS=EISA : XCS=XC : XSS=XS : 
YCS=YC : YSS=YS : ZCS=ZC : ZSS=ZS : FAJDS=FAJDISA : AJDS=AJDISA 
870 GOSUB 940 

880 DELTAISA=DELTAS :R0I=R0S : DISISAIST=DISDS : DISISADST=DISVS : 
ELISTISA=ELEDS : ELDSTISA=ELEVS : DDDI=DDD2 
890 IF ELISTISA >= ANGISA THEN ISEISAIST=1 
900 IF ELDSTISA >= 0 THEN ISEISADST=1 
910 XRI=XR: YRI=YR: ZRI=ZR 

920 ELDSTISA=ELEVS : GMISTISA=GAMAD : GMDSTISA=GAMAV 
930 GOTO 1280 

940 REM SUBROUTINE COM1+COM2 :::::::::: 

950 REM SUBROUTINE COM1 : : : 

960 ETA=SQR (398613.52 #/ AAS ) / AAS : AM=ETA*TMTPS : E=PI 
970 FUN=E-EPSIS*SIN(E) 

980 IF ABS (FUN-AM) <1E-10 THEN GOTO 1030 

990 IF FUN>AM THEN E=E-DEL 

1000 IF FUN<AM THEN E=E+DEL 

1010 DEL=ABS( FUN-AM)/ 2 

1020 GOTO 970 

1030 ROS=AAS*(1-EPSIS*COS(E) ) 

1040 SINPHIO=AAS*SIN(E) *SQR ( 1-EPSIS A 2 ) /ROS : 

COSPHI 0=AAS * (COS (E) -EPSIS)/R0S 
1050 XO=ROS*COSPHIO : Y0=R0S*SINPHI0 : Z0=0 

1060 XI=XCS*X0+XSS*Y0 : YI=YCS*X0+YSS*Y0 : ZI=ZCS*X0+ZSS*Y0 : 

ZZZ= (EPSIS*SINPHIO) / (SQR ( 1+2*EPSIS*COSPHIO+EPSIS A 2 ) ) 

1070 DELTAS=ATN ( ZZZ/SQR ( 1-ZZZ A 2 ) ) : ISIGN=1 : IF FAJDS<.5 
THEN ISIGN=-1 

1080 TC= ( INT (AJDS) + . 5*ISIGN-2415020 ! ) /36525 ! :ALFA=99 . 6909833#+ 
36000 . 7689#*TC+3 . 8708E-04*TC*TC 
1090 OMETE=ALFA+. 25068447#* ( (FAJDS- . 5*ISIGN) *1440!+TIME) : 
OMETE=OMETE/ 360: OMETE= ( OMETE-INT ( OMETE ) ) * 3 6 0 
1100 OMETE=OMETE * CONVERT :XR=XI* COS (OMETE) +YI*SIN (OMETE) : 
YR=-XI*SIN (OMETE) +YI*COS (OMETE) : ZR=ZI 
1110 REM SUBROUTINE COM2::: 
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1120 ZZZ=ZR/SQR(XR A 2+YR A 2+ZR A 2) : SLAT=ATN (ZZZ/SQR ( 1-ZZZ A 2 ) ) 

1130 IF XR=0 THEN GOTO 1190 

1140 IF (XR>=0) AND (YR>=0) THEN SLON=-ATN ( YR/XR) 

1150 IF (XR<=0) AND (YR>=0) THEN SLON=PI+ATN ( YR/ABS (XR) ) 

1160 IF (XR<=0) AND (YR<=0) THEN SLON=PI/2+ATN (ABS (XR/YR) ) 

1170 IF (XR>=0) AND (YR<=0) THEN SLON=ATN (ABS ( YR) /XR) 

1180 GOTO 1210 

1190 IF YR>0 THEN SLON=-PI/2 
1200 IF YR<0 THEN SLON=PI/2 

1210 ZZZ=COS (DLATS) *COS(SLAT) *COS (SLON-DLONS) +SIN (DLATS) * 

SIN (SLAT) : COSGAMAD=ZZZ : GAMAD=PI/2-ATN (ZZZ/SQR (1-ZZZ A 2) ) 

1220 DISDS=ROS*SQR ( 1+ (RE/ROS) A 2-2* (RE/ROS) *COSGAMAD) : ZZZ= (DISDS A 2 
+RE A 2— R0S A 2 ) /2/RE/DISDS : PSID=PI/2— ATN (ZZZ/SQR (1— ZZZ A 2) ) 

1230 ELEDS=PSID-PI/2 

1240 ZZZ=COS (VLATS) *COS(SLAT) *COS (SLON-VLONS) +SIN (VLATS) * 

SIN (SLAT) : COSGAMAV=ZZZ : GAMAV=PI/2-ATN (ZZZ/SQR ( 1-ZZZ A 2 ) ) 

1250 DISVS=R0S*SQR(1+ (RE/ROS) A 2-2* (RE/ROS) *COSGAMAV) :ZZZ1= 

( DI SVS A 2 +RE A 2 -ROS A 2 ) / 2/RE/ DI S VS : PS I V=PI/ 2 -ATN ( Z Z Z 1/ 

SQR (1-ZZZ1 A 2) ) 

1260 ELEVS=PSIV-PI/2:DDD2=SQR(2*RE A 2* (1-ZZZ) ) 

1270 RETURN 

1280 IF (ISEISAIST=1) AND (ELISTISA<ANGISA) THEN ISEISAIST=0 
1290 ISEE=ISEISAIST*ISEISADST*ISEDSADST 

1300 REM INTERFERENCE TO SIGNAL POWER CALCULATIONS:::::::::: 

1310 IF ISEE=1 THEN GOTO 1370 
1320 IF ICOMT=0 THEN TIME3=TIME 
1330 IF ICOMT=0 THEN IMM=IMM+1 

1340 ICOMT=l : TIME4=TIME : TOTTIME56=TOTTIME56+TIME6-TIME5 : TIME5=0 : 
TIME6=0:ICONT=0 

1350 REM PRINT "TIME 3 = "CSNG (TIME3 ) ; "TIME4 = "CSNG (TIME4 ) ; 

"TOTTIME56 = "CSNG (TOTTIME56) 

1360 IF ISEE=0 THEN GOTO 1850 
1370 IF ICONT=0 THEN TIME5=TIME 
1380 IF ICONT=0 THEN INN=INN+1 

1390 ICONT=1:TIME6=TIME:TOTTIME34=TOTTIME34+TIME4-TIME3:TIME3=0: 
TIME4=0:ICOMT=0 

1400 REM PRINT "TIME 5 = "CSNG (TIME5) ; "TIME6 = "CSNG (TIME6) ; 
"TOTTIME34 = "CSNG (TOTTIME34 ) 

1410 ZZZ= (DISISAIST A 2+DISISADST A 2-DISDSTIST A 2 ) /2/DISISAIST/ 
DISISADST:THISADST=PI/2-ATN( ZZZ/SQR ( 1-ZZZ A 2) ) 

1420 DIST3=SQR ( (XR2-XR) A 2+(YR2-YR) A 2+(ZR2-ZR) A 2) 

1430 ZZZ= (DISDSADST A 2+DISISADST A 2-DIST3 A 2 ) /2/DISDSADST/DISISADST: 
THDSTISA=PI/2-ATN (ZZZ/SQR ( 1-ZZZ A 2 ) ) 

1440 REM COMPUTING GAIN OF THE DESIRED GROUND STATION:::::: 

1450 IF THDSTISA=0 THEN GV=1 ELSE GV= (SIN (THDSTISA) /THDSTISA) A 2 
1460 REM COMPUTING GAIN OF THE INTERFERING SATELLITE :::::: 

1470 IF THISADST=0 THEN GD=1 ELSE GD= (SIN (THISADST) /THISADST) A 2 
1480 GV=GNDST*GV:GD=GNISA*GD 

1490 REM COMPUTING DOPPLER EFFECT INCLUDING THE EARTH ROTATION : : 
1500 ZZZ= (DISISADST A 2+ (ROI-RE) A 2-DDDI A 2 ) /2/DISISADST/ (ROI-RE) : 
THT=PI/2-ATN (ZZZ/SQR (1-ZZZ A 2) ) 

1510 DSTLONC=DSTLON-CONVERT* . 1 : ZZZ=COS (DSTLAT) *COS (SLAT) * 

COS (DSTLONC-SLON) +SIN (DSTLAT) *SIN(SLAT) : 

DOP=ROI *SQR ( 1+ ( RE/ROI ) A 2-2*ZZZ* (RE/ROI) ) 

1520 DNP=RE*SQR (2* (1-COS (DSTLAT) *COS (DSTLAT) *COS (CONVERT) - 
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SIN (DSTLAT) *SIN (DSTLAT) ) ) 

1530 COSANG= (DISISADST A 2+DNP A 2-DOP A 2 ) /2/DISISADST/DNP 

1540 VELO=SQR ( 398613 . 52#* (2/R0I-1/AISA) ) *COS (DELTAISA+PI/2-THT) + 

. 463603*COSANG 

1550 FR= FRIST* ( l+VELO/300000 ! ) / (SQR( 1- (VELO/3 00000 ! ) A 2) ) : 
DOPPLER=FRIST-FR 

1560 REM COMPUTATION OF FREQUENCY OFFSET FACTOR FOF 
(INCLUDING DOPPLER EFFECT) :::::: 

1570 DELF=ABS (FRDST-FR) :UPLMT=DELF+BW/2 :UPLMT=PI*UPLMT/AISTBR: 

LRLMT=DELF-BW/2 : LRLMT=LRLMT*PI/AISTBR: DF= (UPLMT-LRLMT) /NDIV 
1580 SU=0 : FOR 1=1 TO (NDIV-1) STEP 2 : F=I*DF+LRLMT 
1590 GOSUB 1670 

1600 SU=SU+4*RESULT: NEXT: FOR 1=2 TO (NDIV-2) STEP 2 : F=I*DF+LRLMT 
1610 GOSUB 1670 

162 0 SU=SU+2 *RESULT : NEXT : F=LRLMT 
1630 GOSUB 1670 
1640 SU=SU+RESULT : F=UPLMT 
1650 GOSUB 1670 

1660 SU=SU+RESULT : FOF=DF*SU/3/PI : GOTO 1710 

1670 REM POWER SPECTRAL DENSITY OF THE INTERFERING MODULATION: : : : 
1680 IF F=0 THEN GOTO 1700 
1690 RESULT= (SIN (F) /F) A 2 : RETURN 
1700 RESULT=1: RETURN 

1710 AISP= (FRDST*DISDSADST/FR/DISISADST) A 2 *GV*GD* FOF/GNDSA/GNDST : 
AISPDB=10*LOG (AISP) /LOG ( 10) :REM PRINT"ISP = "AISPDB ? " (dB) " 
1720 SUMDB=SUMDB+AISPDB: IDB=IDB+1: IF IDB>1 THEN GOTO 1740 
1730 AIMAX=AISPDB:AIMIN=AISPDB 
1740 IF AISPDB<AISPLL THEN GOTO 1770 
1750 IF AISPDB>AISPUL THEN GOTO 1780 
1760 GOTO 1790 

1770 IPDF (0) =IPDF (0) +1 : IJK=0 :MSPECIAL=MSPECIAL+1 : GOTO 1800 
1780 IPDF ( 1001) =IPDF ( 100) +1 : IJK=1001 : NSPECIAL=NSPECIAL+1 : 

GOTO 1800 

1790 IJK=INT( (AISPDB-AISPLL) /DELAISP) : IPDF ( IJK+1) =IPDF ( IJK+1) +1 
1800 MCOUNT=MCOUNT+ 1 

1810 IF AISPDB>AIMAX THEN AIMAX=AISPDB 
1820 IF AISPDB<AIMIN THEN AIMIN=AISPDB 

1830 TMTPISA=TMTPISA+DELT/10 : TMTPDSA=TMTPDSA+DELT/10 : TIME=TIME+ 
DELT/ 10/60 :TIME2=TIME2+DELT/ 10/60 
1840 GOTO 730 

1850 TMTPISA=TMTPISA+DELT:TMTPDSA=TMTPDSA+DELT:TIME=TIME+DELT/60: 
TIME2=TIME2+DELT/60 : IPDF ( 0) =IPDF ( 0) +10 
1860 GOTO 730 

1870 REM COMPUTING STATISTICAL QUANTITIES ::::::: 

1880 PRINT: PRINT: PRINT" % POINTS ABOVE THE SPECIFIED ISP UPPER 
LIMIT OF "AISPUL;" (dB) = "CSNG ( 100*NSPECIAL/MCOUNT) 

1890 PRINT" % POINTS BELOW THE SPECIFIED ISP LOWER LIMIT OF 
"AISPLL;" (dB) = "CSNG ( 100*MSPECIAL/MCOUNT) 

1900 PRINT :PRINT"MINIMUM VALUE OF ISP = "CSNG (AIMIN) ; " (dB) " : 

PRINT "MAXIMUM VALUE OF ISP = "CSNG (AIMAX) ; " (dB) ": PRINT 
1910 SUM=0 : SQ=0 : ICOUNT=0 

1920 FOR 1=1 TO 1000 : ICOUNT=ICOUNT+IPDF (I) : DJMID (I) =AISPLL+ 
DELAISP* (1+2*1 )/2:SUM=SUM+IPDF( I) *DJMID(I) : NEXT 
1930 AMEAN=SUMDB/IDB : FOR 1=1 TO 1000 : SQ=SQ+IPDF ( I) * (DJMID ( I ) - 
AMEAN) A 2 : NEXT : VARIANCE=SQ/ ( I COUNT- 1) 
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1940 DE VI AT=SQR (VARIANCE) 

1950 PRINT "MEAN VALUE OF ISP = "CSNG (AMEAN) ; " (dB) " : 

PRINT "STANDARD DEVIATION OF ISP = "CSNG (DEVI AT) ?" (dB) " 
1960 REM COMPUTING PROBABILITY DISTRIBUTION OF ISP: : : : : : 

1970 SUM=0 : IF ICOUNT=0 THEN GOTO 2000 

1980 FOR 1=1 TO 1000 : POFY ( I ) =IPDF ( I ) /ICOUNT : 

PDF ( I ) =POFY ( I ) /DELAISP : NEXT 

1990 CDF ( 1 ) =POFY ( 1 ) : FOR 1=2 TO 1000 : CDF (I) =CDF (1-1) +POFY (I) : 
NEXT: GOTO 2010 

2000 PRINT "THERE WAS NO INTERFERENCE :::::":GOTO 2070 
2010 TOTTIME56=TOTTIME56+TIME6-TIME5: 
TOTTIME34=TOTTIME34+TIME4-TIME3 
2020 PRINT : PRINT "AVERAGE DURATION OF INTERFERENCE = 

"CSNG (TOTTIME56/INN) ; "MINUTES" 

2030 PRINT "AVERAGE DURATION FREE OF INTERFERENCE = 

"CSNG (TOTTIME34/IMM) ; "MINUTES" : PRINT 
2040 INPUT ABC 

2050 PRINT : FOR 1=0 TO 1001:PRINT CSNG(DJMID(I) ) ;CSNG (IPDF ( I) ) 
CSNG (POFY (I) ) ; CSNG (PDF (I ) ) ;CSNG(CDF(I) ) :NEXT 
2060 OPEN "O" , # 1 , FILEA$ : FOR 1=0 TO 1001:PRINT #1,DJMID(I), 
IPDF ( I ) , POFY ( I ) , PDF (I) , CDF (I) :NEXT 
2070 END 
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